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Abstract 

We introduce and analyse a class of fragmentation-coalescence processes defined on finite systems 
of particles organised into clusters. Coalescent events merge multiple clusters simultaneously to form 
a single larger cluster, while fragmentation breaks up a cluster into a collection of singletons. Under 
mild conditions on the coalescence rates, we show that the distribution of cluster sizes becomes non- 
random in the thermodynamic limit. Moreover, we discover that in the limit of small fragmentation 
rate these processes exhibit a universal cluster size distribution regardless of the details of the rates, 
following a power law with exponent 3/2. 


1 Introduction 


fragmentation, have been widely studied in physical chem- 
19 a century ago (see e.g. 1,12 


Processes of coalescence and its reverse, 
istry since the seminal work of Smoluchowski 
therein). Aside from chemical systems, these processes also have an important role to play in mod- 


18 and references 


elling genealogy 13 20j and even the dynamics of group formation [8j. The microscopic specification of 
a fragmentation-coalescence process is usually given in terms of a stochastic process acting on a finite 
number of constituent particles, introduced by Marcus 16 , Lushnikov [14 15 and Gueron |Yj. However, 


in applied work is it common to write deterministic ‘mean-field’ equations that are intended to describe 
behaviour in thermodynamic limit of large system size. An important question arises: are these processes 
self-averaging so that mean-field calculations are relevant |6,11 


Another key concept in the understanding of large-scale interacting systems is that of universality — 
that certain important macroscopic properties often do not depend on the detailed features of the particles 
and dynamics involved, but rather a much smaller set of properties determine how these processes behave 
in the thermodynamic limit. This is particularly true of critical phenomena, a famous example being 
the directed percolation phase transistion 10 . In fragmentation-coalescence processes, the antagonistic 


nature of the driving mechanisms can give rise to self-organised criticality in certain limits, whereby the 
system evolves to a stationary distribution that exhibits scale-invariant behaviour normally characteristic 
of a phase transition 


17 . It is therefore natural to ask if this behaviour is universal. 


Here we answer both of the above questions in the affirmative. We examine a class of fragmentation- 
coalescence processes that allow simultaneous coagulation of multiple clusters with different rates, com¬ 
bined with independent fragmentation events that shatter clusters to singletons. The fragmentation and 
coalescence rates are taken to be independent of cluster size, and we make a mild technical assumption 
on the coalescence rates. We show that processes of this class approach a non-random limit as the system 
size grows, and moreover, that in the limit of small fragmentation rate we observe a universal cluster 
size distribution with a power-law tail of exponent 3/2. 
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1.1 Finite Fragmentation-Coalescence Processes 

Consider a collection of n identical particles, grouped together into some number of clusters. We define 
a stochastic dynamical process as follows: 

1. Every subset of k clusters coalesces at rate a(k)n 1 ~ k , independently of everything else that happens 
in the system. The coalescing clusters are merged to form a single cluster with size equal to the 
sum of the sizes of the merged clusters. 

2. Clusters fragment at rate A„ > 0, independently of everything else that happens in the system. 
Fragmentation of a cluster of size k results in k ‘singleton’ clusters of size one. 

The standard choice of initial condition is the state with n singleton clusters and this will be the case 
throughout this paper. The factor of n l ~ k appearing in the coalescence rates is included to compensate 
for the combinatorial explosion in the number of subsets of k clusters as n gets larger. In this way, when 
there are order n clusters, the global rates of fragmentation and coalescence of any number of clusters 
are all order n (in the case where X n = A £ (0,oo)). Note that this choice is necessary to ensure that 
there is a single dominant timescale for the dynamics. 

For each n £ N, and k £ {1,..., n}, the state of the system is specified by the number of clusters of 
size k at time t. To that end we introduce the random variables 

w n k{t ) := — #{clusters of size k at time t}, 1 < k < n, 

n 

which we take to be continuous from the left, with right limits, and combine these into the random vector 

w n (t) := (w„,i(t),..., »„,„(!)), 


which takes values on the simplex 


Sr, — 


w r 


= (i 


? n,I? ■ 


e {i, 


i} for i £ { 1 , 




fc=l 


kw. 


n,k 


= l 


Another natural quantity is the empirical cluster size distribution, defined by 

#{clusters of size k at time t} 

-#{clusters at time (}-' 1 “ * “ " 


(1) 


From the description of the dynamics at the beginning of this subsection we can write down the infinites¬ 
imal generator of w n (t ) by summing over all possible events. To do this we must consider the number 
of possible ways of coalescing clusters of sizes l \,..., Ik- If the li are distinct then we simply multiply 
together the number of clusters of each size. If some of the li are the same then there is a complication 
as a cluster cannot coalesce with itself. To each configuration li,... ,1^ we can associate a partition ir 
of {1,..., k} by l u = l v if and only if u,v £ 7q for some i. In this way the set Vk, of all partitions of 
enumerates all the ways we could have chosen k cluster sizes with replacement. Hence, 


•Al nf(lVn) A„n^ ^ W n T 1^ f (l^Vn ) ^ ^ n,k 

n /1 \ n n 

+ EE'"E 7n(/,7T,W„), 


( 2 ) 


k—2 


TTdz'Pk ll=l l \ 7r | —1 


where 


|>! 1 \ \ M 

bn(/> 7r 7 w n) := | / ( w n ~ ^ T Hi 1 / i w n) 1 ) | ffi |, 


(3) 


where e n ^ is a vector of length n with one in position i and zero everywhere else, |7r| is the number of 
blocks in the partition, 7r,; is the number of elements in the i th block of n and (y) z is the falling factorial 
Pochhammer symbol. 
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Rather than working with w n directly, considerable simplification is possible using the empirical 
generating functions G n : [0,1] x R + —>• R 


G„(x, t) = x k w n , k {t), n > 1, 

fc =1 


(4) 


and g n : [0,1] x R+ —> R 


n 

9n{x,f) := y^ i x k p rli k(t) 
k =1 


G n (x,t) 

G„(l,i)’ 


n > 1 . 


(5) 


1.2 Main Results 

Our first main result shows that processes in the class outlined in o are self-averaging with respect to 
the distribution of cluster sizes as n —> oo. 

Theorem 1. Suppose that the coalescence rates a : N —> R + satisfy 


a(k) < Cexp(jk\n\n(k)), \/k, ( 6 ) 

for some constants C > 0 and 7 < 1. Define A = lim-n^oo X n , and let G : [0,1] x R + —>■ R be the unique 
solution of the deterministic initial value problem 


dt 


G( x, 0) = x, 
k\ 


dG (x, t ) = A(x - G( X, t)) + jr (G(x, t) k - kG( 1 , t) k ~ 3 G{x, t)) . 


(7) 


k—2 


If A > 0 then in the limit n —> oo the empirical generating function G n (x,t) defined in W converges 
to G(x,t) in L 2 , uniformly in x and t. That is, 

sup E [(G(x, t) — G n (x, t)) 2 ] < —j= — > 0, as n —> 00 , 

xe[o,i],teR+ V n 


for some constant H. If A = 0 then the above result holds for t less than any fixed finite time T. 

Remark 1. It should be noted that condition © is a sufficient technical condition, but may not be 
minimal. In addition, it is not known if n -1 ^ 2 is the optimal rate of convergence. 

Unpacking the limiting generating function G{x, t) allows us to write Smoluchowski-type differential 
equations for the (rescaled) number of clusters of a given size. 

Corollary 1. Under the assumptions of theorem^ we have L 2 convergence w n j(t) —> Wj(t), obeying 


w j(°) = S j, i) 

f t wi{t) = A(1 - wi(f)) - W X {t) 

k—2 ^ ' »-i 


k-1 


<Z=1 


k-1 


( 8 ) 


d (t) = - a Wj (t) - Wj (t) ( T,Mt) j w iAt)---w lk {t). 

=9 ^ '* \l =1 / k—2 ' h 

h~\ - \-lk=j 


— Wj 

dt 3 


k—2 


for j > 2. Furthermore, the limiting cluster size distribution may be obtained viapkft ) = Wk(t)/ w i(t). 

Our second main result deals with how the limiting process behaves in long times when fragmentation 
events are rare. To that end, we define the stationary cluster size distribution 


p k := lim lim p n ,k(t), 

£—>•00 n —>00 

where p n ,k{t ) is defined in (jTJ) and convergence in the large n limit should be understood in terms of 
Theorem [lj 
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Figure 1: Example of a fragmentation-coalescence process with fragmentation rate A = 10 ~ 7 and co¬ 
alescence rates a(3) = 1, a(4) = 2 and a(k) = 0 otherwise. Left: illustration of possible events (a) 
coalescence of three clusters, (b) coalescence of four clusters, (c) fragmentation of a cluster. Right: com¬ 
parison between a single simulation run with n = 10 9 particles, stopped at time T = 10 6 (black) and 
the n —> oo limit given by Theorem [2] (pale blue). The higher points correspond to the odd sized blocks, 
the lower ones correspond to even sized blocks which are asymptotically vanishing in the limit n —» oo, 
A 0. 


Theorem 2. If a satisfies ([b]) and m is the smallest integer such that a(m) > 0, then the stationary 
cluster size distribution exists and is unique. In particular, the limit of this stationary distribution as 
A —» 0 exists and obeys 


lim p k 

x\o 



if m — 1 divides k — 1 
otherwise. 


In particular lim^o Pk ~ k 3 / 2 , if m — 1 divides k — 1, where f ~ g means f{k)/g{k) — ► c for some 
positive constant c, as k —>■ oo, regardless of a. 

The power-law tail of cluster size distribution is reminiscent of critical behaviour occurring at the 
gellation time in some pure coalescent processes [l], although here it is not a transient phenomenon 
but rather is approached in the long time limit. This is the essential characteristic of self-organised 
criticality, however, there is an unusual complication here that the result depends on the strict order of 
the n — > oo and A — > 0 limits; taking A„ -> 0 we recover a pure coalescent process that does not posses 
this characteristic power-law. 

While it is quite typical for critical processes to exhibit universality in the scaling exponents, our 
result says something stronger that in fact the limiting cluster size distribution is almost completely 
independent of the coagulation rates, depending only on the identity of the first non-zero rate. A 
surprising consequence is that if, for example, the model coalesces clusters in groups of three and four 
(but not pairs) then in the large n and small A limit we will see no clusters of even size whatsoever in 
the stationary distribution. Figure [T] shows an example of this phenomenon for the model with rates 
a(k) = 5k, 3 + 15k, a ■ The model has the apparently paradoxical feature that clusters of even size are 
vanishingly rare, despite the fact that some two-thirds of clusters are singletons, and a(4) > a(3). 

We give here a brief heuristic as to why this is the case. The coalescence rates are scaled so that 
when the number of clusters is of order n the total coalescence rate of events involving k blocks is of the 
same order for all k. However, Theorem [2] states as A —► 0 that the stationary cluster size distribution 
follows a power law with exponent 3/2, which suggests that the average cluster has size of order n 1 / 2 . 
Hence, the number of clusters falls below order n. This causes coalescence events of the smallest number 
of clusters, m, to occur much more frequently than all others and so eventually we only see clusters that 
could have arised from coalescing m clusters. These are clusters of size k where m — 1 divides k — 1. In 
the above example where m — 3, this results in the distribution only giving mass to odd-sized blocks. 
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1.3 Examples of related work 

Without fragmentation one is left with a pure coalescent process 1 . The prototypical example of this 
class has pairs of clusters coalescing at unit rate. In this case, the w n , k (t) converge in probability to 
solutions to the Smoluchowski coagulation equations 19 with unit rate kernel: 


d 1 ^ 00 

= 2 H - w k (t ) 

j =1 3= 1 


Wj(t), 


k > 1. 


These equations can be solved to show 

Wk{t) = 




-2 


2 + t 


k -1 


k> 1, f>0, 


and so the proportion of clusters of size k decays exponentially in k. This behaviour is typical of pure 
coalescence processes. 

A simple modification of these dynamics is to sized-biased fragmentation-coalescence processes, a 
particular example of which (again with only pairs coalescing) has been studied in the context of terrorist 
networks [I]. A thorough analysis of this model was undertaken by Rath and Toth in 17 . For a constant 
fragmentation rate A it was found that the w ntk (t), converge in probability as n —> oo to w k (t) solving a 
set of differential equations similar to those of Smoluchowski: 


1 


k-1 


-^Wfc(t) = ~ j)wj(t)w k -j(t) - kw k (t) - A kw k (t), 

1 3 =1 


k > 2, 


j t w i(t) = —wi(t) + X^^j 2 Wj(t). 

3 =2 

And so it was shown in the long time stationary limit, as t —> oo, that 

A 2 


w k 


exp k log 1 — 


(1 + A) 5 


fc- 5 / 2 . 


As in our model a power-law tail emerges in the limit A —> 0, this time with exponent 5/2. Interestingly 
in this case the critical behaviour survives the joint limit A n —> 0 as n —> oo; the main effort of 17 was 
focused on this scaling regime. 

The critical behaviour of the stationary cluster size distribution, found in Theorem [2j has also been 
seen in a different form of fragmentation-coalescence process studied by Bressaud and Fournier 5|. This 
was a mean-field approximation to a one-dimensional forest-fire model on Z. Edges (each of weight one) 
are glued together if the site between them is occupied, which happens at rate 1, and clusters of k edges 
have all sites removed at rate (k—l)/n. A mean-field approximation to this model, ignoring correlations, 
follows 

4™i n) 0) = ~2w[ n \t) + ^V'fc(fc- 1 )w[ n \t), 


dt 


k> 2 


-1 


k—1 


-<'(i) = -(2 + (fc-l)/nK 


(n) 


W+ 5I w \ n) ^> 55 Wi n) (o«4_i(o> k ^ 2 - 


. i> i 


i=l 


Then, they found that the stationary cluster size distribution converges weakly to (p k ) k >i as n —> oo 
(‘fragmentation’ rate converges to zero) where 

2 (2k-2\ 

Pk ~ 4Nt \ A: — 1 ) 

which matches the distribution in Theorem [2] in the case where m = 2. 

Finally we mention an important separate class of fragmentation-coalescence processes. The models 
discussed so far all concern large but finite particle systems, but it is possible to define (non-size-biased) 
fragmentation and coalescence processes on the partitions of N, as introduced in |fi]. The class of processes 
we work with here can be seen as a natural finite counterpart to these models. 
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2 Thermodynamic limit 


We will focus our proof on the case of fixed fragmentation rate \ n = A £ (0, oo); the extension to the 
joint limit differs only in a few places. Besides standard generating function technology, our methods 
rely on Gronwall’s Inequality [22] which we reproduce below for convenience. 

Lemma 1 (Gronwall’s Inequality). Suppose f : [0,oo) —> ffi. is differentiable and satisfies the following 
differential inequality 

< af(t)+ b, t> 0, 

where a,b £ R. Then, we have that for all t > 0 


m< 


/(0)e Qt + b - (e a 
/(0) + bt 


1) ifa^ 0 

if a = 0. 


The idea of the proof of Theorem |T] is to bound the derivative of the expected squared difference 
between G{x , t) and G n {x , t) in such a way so as we may apply lemma [l] where, here, a will be negative 
and b = b n will decay to zero as n —> oo. To achieve this, we will need to prove that the derivative exists 
and certain error quantities (specified in the next section) converge to zero as n tends to infinity. 


2.1 Mean-Field Calculation 


We first undertake a mean-field calculation to determine a viable candidate for limn-^oo G n under self¬ 
averaging. 

Lemma 2. For f 0 (x,w n ) := ^L-=i xkw n,k> we have 

n a(k) 

A n f 0 (x,W n ) = \(X - / 0 (x,W n )) + Y] fp{x,W n ) k - kfo{l,W n ) k ~ 1 f 0 (x,W n )) + /3 n (x,W n ), 

k—2 


where 

A 

sup \Pn(x,W rl )\ < —, 

where A is a constant independent of n and x. 

Proof. We use the infinitesimal generator from ([2]) in the special case where / = /o using a slight abuse 
of notation. Hence 


Anfo {x, 'Wn) 


n n 

-A ^2 x k vJ n ,k + A ^2 kxw n ,k 

k=1 k =1 


+E 


fc=2 


a(k) 

k\n k 


n n 


EE-E 

ir£Tkh=l *M =1 



kil h 


kl \ kl 

El 7Ti \ xU Yl^ nWn di)\^i\- 
i=l J i=1 


As n grows large, we claim that the dominant contribution from the final product inside the third term 
on the right-hand side is simply • • • w n j k . If n = {1,..., fc}, then this is the only term, otherwise 

there is a subdominant correction resulting from the fact that clusters cannot coagulate with themselves. 
Labelling this correction as P n (x,w n ), which we will later bound, we obtain 


Anfo( x , w n) = “A ^2 xkw n,k + A ^ kx W n ,k 


k=l 


1 

n n 


/ j \ n 'n / k k \ 

TlZ-EP-PK-r*.* 

k= 2 q = l l k = 1 \i=l *=i / 

"h Pn (x, W Tl ) 

*{k) 


= Hx - fo(x,w n )) -jfp- ifo(x,w n ) k - kfo(l,w n ) k 1 fo(x,w Tl ))+p n {x,w ri ). (9) 


k—2 
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It remains to show that f3 n (x,w n ) < A/n for some constant A > 0. We know it is the remainder 
after subtracting the dominant term and so is 


i3 n {x,w n )=Y J y y---y ( n^ 11 *i in 


M 


M 


M / \ M 

{nw n>h ) ki | 


fe =2 


k\ 


TT^Vk 11=1 l |tt|=1 \ i——l 


i—1 


1 


ik,l 


-n»- 

i= 1 


Each multinomial in the second bracket, when expanded, has at most 2 k terms with varying powers of 
1/n larger than or equal to one. Also, we can see that all the terms are divisible by Wn.q • ■ • and 

have a coefficient that is less than or equal to (k — 1) J , where j is the power of 1/n in that term. Hence, 

, a . ,.^-A a(k) v- V- A 2 fc (fc-l) 

| Pn{x,Wn)\ < 2 J ,j, jp / ,, / y 2 ^ n Wn,Ji ' ' ' • 

fe=2 ^ TreP), ii = l Z|„| =1 

Summing over Zi,..., li n i we have that, 

\Pn(x,w n )\ < MW’ 1 - 

fc=2 ' ' 7T Gpfc 


Then noting that /o(l,iUn)^ < 1 for all tt , we have 


lA a(k)2 k „ 

|/?n(£,«V)| < n ^2^_ 2 yBk, 


where Bk denotes the k th Bell number and is the size of the set Vk- We bound Bk using a recent result 
from [ 2 ], obtaining 


|/3 n (x,w n )| < - 


(h - 


a(k) ( 1.584A; 


k=2 


(k — 2 )! \ln(fc + 1 ) 


If we look at the asymptotic behaviour of the summands for large k, we see that 

k 

~ exp (— k In ln(A) + 0(k )), 


1 


1.584& 


{k - 2)! \\n(k+ 1) / 

where ~ is defined in Theorem [ 2 ] Thus, by assumption ([b]) on a, the series converges as required. □ 


2.2 Self-averaging 

With the mean-field behaviour determined, we proceed to the proof of L 2 convergence for sample paths. 
The following two lemmas will help us determine the behaviour of G 

Lemma 3. The function fo(x,w n ) = Y^k=i xkw n,k satisfies 

Xx 2 

A n f() (x, — 2fo(x, W n )Anfo(x, 'fXn) T ~ f ‘2 (1; ^Jn ) T tl(x } 


where 


f 2 (x,w n ) := y k 2 x k w n k, and sup h(x,w n ) < —, 

z r* TL 


k =1 


IV n E«Sn 


such that E is a constant independent of n and x. 

Proof. We plug the function f§ into the infinitesimal generator ([ 2 ]) 


n ^ 1 / kx — x k \ \ 

^nfo (*^5 ^n) — ^ ^ | ( fo(%i ^n) “ 1 “ ~ J fo j Wn,k 


k =2 7TE'Pfc il = l 1\tt\ —1 
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where 7 is defined in ([ 3 ]). By expanding out the brackets we see that 

Xx' 2 

Anfo ^n) — 2/o (X, IV n ) An f 'o {x , IVnfj “b — /^(l; “t“ h[x, W n fi 

where 


h(x, w n ) = - (x 2k - 2kx k+1 ) w n , k 


k =1 


;tf EE-E (IB'"'“-£m*‘ 

fc=2 K£V k l 1= l Zi„i =1 \i=l 


M 


M 


i=l 


M 

n 

2=1 


(nw nth ) | ff .| 


)>d 


< -—a ; 2 + -fo(x 2 ,w n ) + - ^ -j-P-fo(x 2 , w n ) k 
n n n k\ 

k—2 

- - ^>2 ~Tr k fo( x 2 , w n)fo(x,W Tl ) k - 1 + - ^ -TT^ fc /o(a; 2 ,^)/o(l,'W „) fc_1 
n * —' ft! n ' f£\ 


k—2 


k—2 


ly^ a(/c) , 2 , , \k—2 1 a(/c)/c / 1.584fc \ fc 

k—2 k—2 

where the bound was attained by using Lemma [4] and the same techniques as in the proof of lemma [2j 
Now, using the bound fo{xfiw n ) 1 < 1 for all j,l > 1, and (10), we obtain 


1 / n / 1 3 

h(x,w n )<- A + J>(fc) ry + 77 — 


k—2 


k\ (k- 1 )! (fc- 2 )! 


B \ E 

H-< — j 

n n 


where E is independent of n and x, as required. 

Now we need to bound the expected value of the function f 2 (l,w rl (t)). 
Lemma 4. Define the function C n : [0, 00 ) -> R as follows 

C„(t) = E [/ 2 (1, w n (t))}, 
for n > 1 and t > 0. Then we have that 


a 


sup C n (t) < D, 
ra£N 


where D is a constant independent oft. 


Proof. We will use Gronwall’s Inequality. We look at the derivative of C n , using Dynkin’s formula we 
see that 

^ C n (t) = E [A n f 2 ( 1, w n (t ))], 


and hence, using the infinitesimal generator ([ 2 ]), 
d C n {t) = K -A 

k= 1 


dt 


- A^(fc 2 - k)Wn : k(t ) 


n / 7 \ n n 

+EfEE-E 

k—2 7T ^1 — 1 ^ I ?r I—1 

which, using a similar argument to that used when bounding fi n (x,t), gives that 


/ M \ 2 M 

»|li 

1- 

-to — 

c* E 

\t=l / 2 — 1 

2=1 J 
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< A — AG„(t) 

at 


where 


Thus, 


•2E 


•2E 


E wr E•••E (E^ ] w n,iAt)■ ■ ■ w n,i k tt) 


n n 


fc! 

k—2 l i = l /fc = l \i<J 


E 

k—2 


a(k) 2 k (k - 1) 
fc! n 


n n 


E E■" E Ei*m I 

we'Pfc ii=i j it|=i \*<j 


< A — XC n (t) + 2E 


< A — A C n (t) + 

k=2 


e (*WM > k| - 2 

,fe=2 ■ fc=2 ■ n£Vk V 7 


a(fc) 13 


(fc — 2)! n’ 


S = E 


a(k)k f 1.584fc 


k—2 


(fc- 2)! \ln(fc + 1) 


( 10 ) 


^<?„(*) < Di - XC n (t) 
at 

where D\ is some constant not dependent on n or t, so applying Gronwall’s inequality gives the result 
as required. 

□ 


With the help of the above three lemmas we are now ready to prove Theorem [l] 

Proof of Theorem 1. Recall from |7| that we write G(x,£) for the solution of the differential equation 


^ Or, t ) = A(x - G(x, t)) + jr ( G( ®, t) k - kG( 1, t) k ~ 1 G(x, t)) 

fc =2 

G(x, 0) = x. 


The existence and uniqueness of G(x, t) is straightforward to establish via the Picard-Lindehof Theorem, 
as the right hand side is uniformly Lipschitz in G, x and t. Specifically, starting with the case x = 1, we 
have 

dG(1,t) =F 1 (G(l,t)), (11) 


dt 


where 


F 1 (g) = \(l-g) + J2 a jr(l-k)g k . 


k> 2 


fc! 


The derivative of Fi is uniformly bounded on g £ [0,1] since 

pt( n \ _ \ a ( fc ) „k -1 


fc >2 


and thus 


sup |i^(ff)| = A + it-l) i = L i < °°- 


se[o,i] “ (fc - 2)! 


( 12 ) 


(13) 


(14) 


Said another way, F\ is uniformly Lipschitz with constant L \, and thus © has a unique solution. A very 
similar argument using furthermore the well-posedness of (111 applies to the case of general x € (0,1], 
where we have 

= Fx(G(x, t), t ), (15) 
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for 


(16) 


Fx{g, t ) = A {x-g) + J2 ( g k - kG(l, t) k 1 g ) . 

k> 2 

The derivative is again bounded: 

K (<?, t) = -A + g (/- 1 - G(l, t)^ 1 ) , 

and since G(l,f) G [0,1] we obtain 

y{k) 


sup \F' x {g,t)\ = A+ 2^ --— =L<oo. 


k> 2 


(17) 


(18) 


As this constant is independent of x, g and i, global existence and uniqueness follows. 

Moving on now to compare to the finite system generating function G n (x , t), we start by noting that 
in the trivial case x = 0, we have 

G n (0,t) = G(0,t) = 0, Vt> 0. 

Now, for strictly positive x we define the function Y n>x (t) : [0, oo) -> R as follows: 

Yn,x(t) = E[(G(x,t) - G n (x, t)) 2 ] , t > 0. 

We will use Gronwall’s Inequality, so we look at the derivative of Y n x (t) with respect to t, using Dynkin’s 
formula (bearing in mind that pre-limit in n, the system is finite) we see that 


dY n 


dt 


*-(t)=E[A n (G{x,t)-G n (x,t)) i 


= 2 G(x, t)g~ t G{x, t)+E [ A n G 2 (x , t)] - 2G(x, t )E [A n G n (x, t)) 

a 

-2E[G n (s,i)]^G(z,i). 


(19) 


Noting that G n (x,t) = fo(x,w n (t)), we may use lemma [2j lemma [3j and the definition of G in |7| to 
see that 

4jYn >x (t) = 2G(x,t)(x(x - G(x,t)) + ^ ~ kG(l,t) k ~ 1 G(x,t)) 

V U—O • / 


2E 


G n (x, t) [X(x - G n (x, t)) + /3 n (x , w n (t)) 

+ t ( G «( ;c ’ - kG "(!> i) k ~ x G n {x, t)) 


k=2 


-2G{x,t)E 


X(x - G n (x,t )) + p n (x,w n {t)) 

+1 t)k ~ kGn{1, ^ k ~ 1<G ^ x ' *)) 


k=2 


- 2E[G„(x, t)] (A(x - G(x, ty + ttr ^ ^ ~ kG{1 ' Q k ~ lG ( x ’ *)) ) 

' k—2 ' ' 


\x^ 

H- C n (t) + E [h(x, w n (t))} 


Ax 2 


= -2A Y n>x (t) + 2E[(G„(x, t) - G( x, t))/3 n (x , w n (t))] + —C n (t) + E[h(x, «>„(*))] 


a(k) 

+ 2}^ -4ae 


k—2 


k\ 


(G(x, t) - G n (x,t)){G(x,t) k - G n (x,t) k ) 
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- A(G(l,t)*- 1 G(a:,t) - G„(l J t) fc - 1 G„(a:,i))) 

+ 2 a ^E[{G{x,t)-G n {x,t)){G(x,t) k -kG(l,t) k - 1 G{x,t)]. (20) 

fc=n+1 


Note that the right hand side of (201 depends on G n (l,t). For the boundary case x = 1, however, we 
have a closed expression in G n (l,f). The plan is thus to first show the theorem holds for x = 1 and use 
this to complete the proof for general x. Define X n (t) := Y n ^(t). 


Lemma 5. We have 


sup X n (t) < 

t> 0 


Hn 

2A 


where 


H n := - 
n 


2A + D\ + E + Y^ 


k—2 


2 a(k) \ 

V-)' 


In particular, H n — > 0 as n —> oo because of assumption ([6]). 


( 21 ) 


Proof. Substituting x = 1 into (201 we see that 


(f) = -2A X n {t) + ^G n (t)+E[h(l,w n (t))] + 2E[(G n (l,t) - G(l,t))/3 n (l,w n (t))} 
+ 2 E ^r E [(^(i,*) - G n (l,t))(G(l,t) k - G„(l, t) k )] (1 - k) 

k—2 

+ 2 E ^n(G(l,t)-G n (l,t))G(l,t) k (l-k). 

k=n -\-1 


Now, (G(l,f) — G n (l,t))(G{l,t) k — G n (l,t) k ) > 0, which means the second line above is negative. Hence, 
together with lennna [3] and lemma [4j we have 


(t) < —2A X n (t) d—(DA + E + 2A) + y'' j- — E-t < —2A X n (t) + H n , 

n k^n+l — X ) 


dt 


so that applying Gronwall’s inequality gives the result as required. 
We can now use this bound to complete the proof of Theorem |T] 


Proof of Theorem 1 (cont.) Continuing from (20), we have that 


( 22 ) 

□ 


~rjYn,x(t) < —2A Y n>x (t) + 2 Sk 

k —2 


Hr, 


(23) 


where 


S k 


a(k) 


E 


G(x,t) k 1 ( G(x,t) 2 


k\ 

- fcG(l,t) fc_1 (G(x,f) 2 — 


- G(x, t)G n (x, t)) + G n (x, t) k 1 (G n (x, t) 2 - G{x , t)G n {x, t)) 
G(x,t)G n (x,t)) - kG n ( l,f) fc_1 ( G n (x,t) 2 - G(x,t)G n (x,t)) 


We will apply lemma[5]to bound the sum of the Sk- First, it is necessary to bound the sum in terms of 
E[|G(l,f) — G n (l,f)|], and remove any terms involving just G(x,t) and G n {x,t). To do this we create 
terms that contain the positive term ( G(x , t) — G n (x , f)) 2 , so that if we pre-multiply them by something 
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negative we can discard it for an upper bound. For example, looking at the first term above we can do 
the following by adding a zero 

G(x,t) fc_1 E [G{x,tf - G{x,t)G n (x,t)\ = G(z,i) fc_1 E [ G(x,t) 2 - G(x,t)G n {x,t)\ 

- G(x,t) k ~ 1 E [G(x,t)G n (x,t) - G n (x,t) 2 ] 

+ G(x, <) fc_1 E \G(x, t)G n (x, t) - G n (x, t ) 2 ] 

= G{x,t) k ~ 1 E |"(G(x, t) — G n (x,t)) 2 

+ G( x, <) fc_1 E [G(x, t)G n (x, t) - G n (x, t ) 2 ] 

= G(a;,t) fc_1 E \(G(x,t) - G n (x,t)) 2 

+ E [G(x, t) k ~ 2 G n (x, t) (G(x, t) 2 - G(x, t)G n (x, t))] . 

In creating the square term above, we get a similar term to what we started with but with the exponent 
of G(x,t) decreased by one and the exponent of G n (x,t ) increased by one. We repeat this process until 
the exponent of G n (x, t) is k — 1. We then subtract and add back in the same terms but with x replaced 
by 1. Specifically, 

S k = ^yy^E G(x,t) k ~ 1 (G(x,t) - G n (x,t)) 2 + G(x,t) k ~ 2 G n (x,t)(G(x,t ) - G n (x,t) 2 ) 
k\ 

+ ... + G{x,t)G n {x 1 t) k ~ 2 {G(x,t) - G n (x, t)) 2 + G„(x,t) k ~ 1 (G(x,t) - G n (x,t)) 2 

- [g(1, t) k -\G(x, t) - G n (x, t)) 2 + G(l, t) fe_2 G„(l, t)(G(x, t) - G n (x, t) 2 ) 
k\ 

+ ... + G(l, t)G n (l, t) k ~ 2 (G(x, t) - G n (x, t)) 2 + G„(l, t) k ~ 1 (G(x, t) - G n (x, t)) 2 
+ ^E G(l,t) fc - 1 (G(*,t) - G n (x, t)) 2 + G(l, t) k ~ 2 G n {l, t)(G(x, t) - G n {x,t)f 
+ ... + G(l, t)G n (l, t) k ~ 2 (G(x, t ) - G n (x, t)) 2 + G„( 1, t) k ~\G(x, t) - G n (x, t)) 2 

- ^ 7 ? (kG{ 1, t) fc - 1 G(®, t) 2 - fcG(l, t) k ~ 1 G(x, t)E[G n (x, t)] 

k! \ 

+ fcE [G n (l, t) fc - 1 G n (*, t) 2 ] - kG(x, t)E [G n (l, tf-'G^x, i)] ) , 

noting that the first four lines combined give something negative. As G(x,t) < G(l,t) and G n (x,t ) < 
G n (l,t), we have 

S k < ^pE [G(l,i) fc " 1 (G(®,t) 2 - 2G(x,t)G n (x,t) + G n (x,t) 2 )\ 

+ ^E [G(l, t) fc_2 G„(l, t)(G(x, t) 2 - 2 G(x, t)G n (x, t ) + G n (x, t) 2 )] 

+ ... + ^E [G(l, t)G„(l, t) k ~ 2 (G(x, t) 2 - 2 G(x, t)G n (x, t) + G n (x, t) 2 )] 

+ ^E [G„(l, t) k ~\G{x, t) 2 - 2G(x, t)G n (x, t ) + G n (x, t) 2 )] 
kG(l,t) k ~*G(x,t) 2 — kG(l,t) k ~ 1 G(x,t)E[G n (x,t)} 

k\ \ 

+ kE [G n (l, t) kJl G n (x, t) 2 ] - kG(x, t )E [G„(l, t) k ~ l G n (x, t)] ^ . 

We gather terms in G(x,t) 2 , G n (x,t) 2 and G(x,t)G n (x,t), each multiplied by \G(l,ty — G„(l,t) J | for 
some j. Here we have taken absolute values in order to not worry about signs. We do this by matching 
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each of the terms from the first four lines of the above with the corresponding term from the last two 
lines. Doing this we see that 

S k < ^G(x,t) 2 E\G(l,t) k ~ 2 |G(M) - G„(M)| + ... + iGfl.t)*- 1 - G n (l,t) fc - 1 | 


jfe! 

+ ^E 

k\ 

a(k) 


G n (x, t) 2 ( |G(l,f) fe_1 - G„(l, t) k ~ l | + ... + G„(l, t) k ~ 2 |G(1, t) - G n ( 1, 


, j^G(x, t)E [c n {x, t) (g(1, t) k ~ 2 |G(1, t) - G„(l, t)\ + ... 

+ G(M)|G(M) fe - 2 -G„(M) fe - 2 |)' 

+ ^rG(x, t)E [G n (x, t) (G n ( 1, t) |G(l, t) k ~ 2 - G n ( 1, t) fe - 2 1 + ... 

+ G n (l,t) fc - 2 |G(l,t)-G n (l,t)|) . 

Hence, using the fact that G(x,t), G(l,t), G n (x,t ), G„(l,t) < 1, we see that 

|G(l,£) fc_1 - G n (l,£) fc - 1 | + ... + |G(1, t) - G n (l,t)\ 


S t <^& E 


jfe! 


*(&), 


+ -^E [(|G(l,£) fc - 1 - G„(l,i) fe_1 | + ... + |G(1, t) - G„(l, f)Q] 

H—^pE [|G(1, £) fc 2 — G n (l,t) k 2 | + ... + |G(l,t) — G„(l,i)|] 

H—jpE [|G(1, 2 — G n (l,t) k 2 | + ... + |G(l,f) — G n (l,£)|] 

= 2^ xiE [|G(1, ty - g„( i, ty |] + 2 ^ x: e O^ 1 - ty - *) J ’|] 


s k < 2 


jfe! 

a(fc) 

~fcT 


j=i 

fc-i 


E-?'E 00(1, t) - G„(l,£)|] + 2°^ Y, JE [|G(l,t) - G„(l, t)|], 


3 = 1 
k—2 


k\ 


3 =1 J=1 

where the above inequality uses Ij/- 7 — | < j|y — z\ for y,z < 1. Thus 

2 a(k){k — 1) 


S k < 


k(k - 2)! 


■E [|G(l,t) — G n (l,f)|] 


We can plug this result back into (23) and see that 

d 


dt 


Yn,x(t ) < -2Ar„,x(t) + 4E [|G(1, t) - G n (l, £)]] ]T Q ' (fc , )(fc 


fc= 2 


jfe (jfe - 2)! 


Now using lemma [5] we have that 


JW n " a( fc)(fc-l) 


dt Yn ’ x{t) ~ 2XYn ’*^ + 4 ^ E _ 2)!' + Hn ' 


Then, we can use Gronwall’s Inequality to see that 

\y/H^^a(k)(k- 1 ) 


Fn,x(£) < 2 


E 


Av^ ^ k(k-2)\ 


+ H n ) (1 - e- 2At ) 


^ 0 V^n p a(fc)(fc- 1) | u 

~ W\h k(k- 2)! + n 

-S> 0 


as n —> 00 , thanks to Q. This convergence is uniform in both x and t, as required. 


(24) 


□ 
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From Theorem 1 the Smoluchowski-type equations stated in Corollary 1 follow immendiately: 
Proof of Corollary 1. We begin by subsituting G(x,t) = Yhj x ^ w j{t) into 

^ (x,t) = A(x — G(x,t )) + ^ (G(x, t) k — kG(x, t)G(l, f) fe_1 ) 


dt 


k =2 


k-l'y 


= X ( x-'EJwM +Y^ 

2=1 / k =2 


i(fc) 

“fcT 


“ M XX'W 


l 2=1 


u =1 


U = 1 


Rearranging the right hand side to collect powers of x we obtain 

^ on 


fc-lN 


9G, N 

~gfM = X 


X{1 - in{t)) - Wl {t)J2 JT^—rX ( X W jW 


V 




a(/c) 


fe-i 




+XX - m- w - w X_ 1V . z_ 

2=2 V fc=l ' ’’ \! = 1 / 

+X^r X w ii(*)■■■«'«*w). 

k=2 / 

h-\ - \~lk = j 

as required, where the final term on the right-hand side is obtained using the product formula for power 
series. □ 


3 Universal Stationary Distribution 


Now that we know how the empirical generating function G n (x,t) converges in the large n limit, we can 
extract a prediction for the stationary distribution of cluster sizes. Recalling from <[5j) that the generating 
function of the cluster size distribution is obtained via 


g n {x,t) = 


G„(x,t) 
G n (l,t) ’ 


we reintroduce (and later prove the existence of) the quantities 


9x(x,t) = 


G(x , t) 

G(ur 


9\{ x ) = lim 9\(x,t) = ’ 

i-> oo Ga(1) 


G\(x) = lim G(x,t ), 

£->-oo 


Pk{ X) = lim lim p n ,k(t )■ 

£—>•00 n—>oo 


Note: as this section focuses on the limit as A —> 0 of the above quantities, these quantities have been 
labelled with A to make it easier to follow. Under self-averaging, the stationary cluster size distribution 
is found simply by determining the fixed point G\(x) of the differential equation ([T]) and developing 
g\{x) = G\(x)/G\(l) as a power series. In practice, this procedure is only tractable in the limit A \ 0, 
where we find a particular limiting distribution. 
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Proof of the convergence of G as t —> oo is an application of the Poincare-Bendixson Theorem (see 
p. 223), which states that attractors of plane autonomous systems are either fixed points, collections of 
fixed points (with their associated linking orbits), or limit cycles. To compute the asymptotic form of 
the solution in small A we employ the classical technique of series inversion that dates back to Newton, 
Lagrange, Biirmann and Puiseux.In particular, the following theorem is a sub-case of results of those 
authors (see [9 , p. 183): 
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Theorem 3 (Series Inversion). Let f : [0, oo) —> [0, oo) be defined by the power series /( x) = Yhk>m a kX k 
for integer m > 1 and non-negative real constants {a*;}. The inverse is expressed via the Puiseux series 

OO 

f~\y) = Y J b k y k/m , 
fc=l 


where 


bk = Jr lim 
k\ z->o 


jk— 1 


dz k ~ l yfiz) 1 / 


In particular, note that b\ = 1/ai. 


Proof of Theorem, 2. Recall that the limit generating function G(x,t) obeys Q. Equivalently we may 
write dG(x,t)/dt = F x (G(x,t),t), where F x is the function defined in (16). Since G(l,t) > G(x,t) for all 
x , it follows that F x is monotonically decreasing in its first argument, moreover, F x (l,t) < 0 < F x (0,t). 
Hence, in the particular case x = 1, G(l,f) has a unique fixed point , denoted G A (1). As G(l,f) £ [0,1] 
for all t convergence is assured. 

For general x, note that (G(l, t), G(x, t )) is a plane autonomous system and, by the same argument as 
above, has a unique fixed point. From the Poincare-Bendixson theorem, we conclude that we need only 
to rule out the case of periodic orbits to confirm convergence to the fixed point. But we have already 
shown that G(l, t) converges to its unique fixed point G A (1), so the w-limit set of the plane system must 
be contained in the line (G A (1), ■)■ Thus, by the Jordan Curve theorem (see 21 p. 220), the w-limit set 
cannot be a regular periodic orbit. Therefore, we may conclude that the solution of ([7]) is required to 
converge to a stationary value G(x,t) —> G x (x), where 


0 = A(x — G\(x)) + jr ^ ( G x (x) k - kG x (l) k ~ 1 G x (x)) . 


k—2 


At the point x = 1 we rearrange to find A as a power series in Ga(1): 

A = rrgyij E ^ir (t " 1)&(1, ‘ = ^ (ra “ 1)Ga(1) ” + E «&(!)' 

k"^ TTL k^>7Tl 


(25) 


(26) 


where m is the least integer such that a(m) > 0, and Ck are the appropriate constants. From Theorem 
[3] we thus determine that 


G a (1) = A 1 /" 1 ( 


TO! 


(to — 1 )a(m) 


l/r, 


+ o( A 1/m ). 


For general x , equation (25) implies the same leading order behaviour of G A (1) in A. Recalling the 
definition g\{x) = G\(x)/G\(l) at the start of this section, we see that 

0 = A (* - G x (l)g x (x)) + £ ^G A (l) fe ( g x (x) k - kg x {x)) , 


= \x — 


1 ! 




(to — 1 )a(m) 
a(k) 


k—m 

1/m 


k—m 


k\ 


X 1+1/m gx(x) 

1 \ k/m 


(to — l)a(m) 


\ k/m ( g x (x) k - kg x (x)) + o(X ), 


= Ax + 


A 


-mg x {x)) + o(X), xG[0,1]. 


TO — 1 

In other words, when to is the smallest integer such that a is non-zero, we find that in the limit as A 0, 

mg x (x) - g x (x) m - x(m - 1) = 0. 
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Note that this limit exists as g\ is a Puiseux series in A. Applying the inversion theorem once more, we 
obtain the series 


OO b, 

x K 


^(*) = E*r~o 


k=1 


jk—1 


dz k ~ x \f{z) 


(27) 


where f(x) = m./(m — 1 )z — l/(m — 1 )z m . We compute 


/(*) 


m — 1 
m 

m — 1 
m 

m — 1 
m 


1 


1 — z m ~ 1 /m 

k / oo 


\n—0 
/c oo 


E - 

eQ 


1 - 1 ) 


n—0 


k — 1 + n 
n 


n{m— 1 ) 


using the product formula for series. To get the summands for (27) we need to look at the (k — 1) 


th 


derivative of the above series, evaluated at z = 0. Hence we have that 
k 

z—0 I 0 


jk—1 


dz k 1 \f{z) / 


\{ w z^) k {b) m 1 {k- !)! if m - l divides k - 1 

\f( z )J 2 =o 1 n otherwise. 


Plugging this into (27) and comparing coefficients of x k , we have 


lim pfc(A) = 
A \0 


1 r m-i \ K / M 
k V m ) V m / 


fc-i 

a-i 


if m — 1 divides k — 1 

^0 otherwise. 

Hence, we can conclude that, for large values of k , under these conditions, 


lim p k {X) ~ k 3/2 , 


regardless of the values of a(k), k> 2. 


□ 
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